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A nanoflare model of quiet Sun EUV emission 
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Abstract. Nanofiares have been proposed as the main source of heating of the solar corona. However, detecting them directly 
has so far proved elusive, and extrapolating to them from the properties of larger brightenings gives unreliable estimates of 
the power-law exponent a characterising their distribution. Here we take the approach of statistically modelling light curves 
representative of the quiet Sun as seen in EUV radiation. The basic assumption is that all quiet-Sun EUV emission is due 
to micro- and nanofiares, whose radiative energies display a power-law distribution. Radiance values in the quiet Sun follow a 
lognormal distribution. This is irrespective of whether the distribution is made over a spatial scan or over a time series. We show 
that these distributions can be reproduced by our simple model. By simultaneously fitting the radiance distribution function and 
the power spectrum obtained from the light curves emitted by transition region and coronal lines the power-law distribution of 
micro- and nanoflare brightenings is constrained. A good statistical match to the measurements is obtained for a steep power- 
law distribution of nanoflare energies, with power-law exponent a > 2. This is consistent with the dominant heat input to the 
corona being provided by nanofiares, i.e., by events with energies around 10 23 erg. In order to reproduce the observed SUMER 
time series approximately 10 3 to 10 4 nanofiares are needed per second throughout the atmosphere of the quiet Sun (assuming 
the nanofiares to cover an average area of 10 13 m 2 ). 
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1. Introduction 

The solar emission in the EUV and X-ray wavelength range 
features transient events on all scales, such as flares and 
micro- or nanofiares, explosive eve nts , and blinkers, see, e.g., 
iTandberg-Hanssen & Emslid (1 1988b or lSolankil (12002b and ref- 
erences therein. Flares, micro- and nanofiares refer to similar 
phenomena with different magnitudes of energy release, which 
are usually studied in coronal emission, while explosive events 
and blinkers are spectral line broadenings and brightenings that 
are mainly restricted to the transition region. Events in the en- 
ergy range 10 30 - 10 33 erg (10 23 - 10 26 J) ar e usually refe rred 
to as "normal" flares. Micro- and nanofiares dParkerll 1988b are 
the brightenings with energy below approximately 10 27 erg, al- 
though the limits vary somewhat in the literature. Such tran- 
sient brightenings are thought to be the signature of dissipation 
of the energy stored in the magnetic field, e.g., through mag- 
netic reconnectio n. Thus, they provide a mechanism for heating 
the solar corona dParkerll 1988b . 

The frequency or rate distribution of the energy of flares 
(dN/dE) has been found to obey a power law for several 
wavelength regimes. Powe r laws have also been found for 
blinkers by Brko vic et al.l (200ol) and for explosive events 
by IWinebarger et al.l (12002b . The distributions found by 



Winebar ger et al.l (|2002), however, look more consistent with 
a lognormal function than a power law. Flares and their dis- 
tribu tions have been the s u bject of a large number of studies, 
e.g.. iDatlowe et all dl974: iLin et al.l dl984h: iKrucker & Bend 
1998b: IParnell & JupdI J2OO0I): lAschwanden et al.1 (|2000h : 
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Mitra-Kraev & Benz (2001); Aschwanden & Parnell (2002); 



Vekstein & Jainl (12003b . Power laws have also been applied fo r 
stellar flare energy d i stributions, see, e .g., Audar d et al. (2000); 



Gtidel et al.l d2003b: IXrzner & Gtidell d2004 and references 
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therein. Lu & Hamilton! ( 1 99 ll) have explained the power-law 
dependence of the solar flare occurrence rate in a model of self- 
organized criticality as avalanches of many small reconnection 
events. 

The observed exponents a in the power-law relation dN/dE = 
E{)E~ a range from 1.5 to 2.9 for solar and for stellar flares, but 
most studies give values of a < 2. Exponents within this range 
are also obtained for blinkers and explosive events. The wide 
range of the exponent is to a large part due to different selec- 
tion criteria or assumptions about geometry and structure of the 
emitting plasma as well as instrumental resolution. The larger 
the exponent, the more weight is given to small-scale events 
such as micro- and nanofiares. For an exponent greater than 2, 
the energy content is dominated by the small-scale events, and 
in order to have finite total energy content, a lower cutoff in 
energy has to be introduced. Since the energy released in flares 
and flare-like events of sufficient size to be directly detected is 
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not enough to explain the temperature of the corona, an expo- 
nent above 2 is necessary if a mechanism producing flaring is 
to be the main cause of coronal heating. Thus, the distribution 
of the smallest flares is of great interest in order to have reliable 
estimates on their participation in coronal heating. 
Most previous work has concentrated on identifying individual 
events in radiance time series. Clearly, it is easy to miss or un- 
derestimate the number of the weakest such events, which may 
lead to a background of almost continuous brightening if they 
are sufficiently common. Note that in order to judge the amount 
of energy released through transient events, the X-ray, EUV 
and UV emission is often decomposed into a (nearly) steady 
background of emission with superposed transient brighten- 
ings. 

In this work we start from the basic assumption that all the 
emission at EUV wavelengths is produced by transient events 
(i.e., the observed apparent background is a superposition of 
many such events). We model in a simple way time series of 
the radiance for a set of parameters that is consistent with ob- 
servations. The results are then statistically compared to the 
observational data. 
The main aims of this paper are 

I. to identify diagnostics which allow such a statistical com- 
parison to be made. The two diagnostics found and studied are 
the lognormal probability density function (PDF) of EUV ra- 
diances and the wavelet or Fourier power spectrum of radiance 
time series; 

II. to test to what extent the use of these diagnostics can con- 
strain the parameters of the model, in particular the power-law 
exponent a of the distribution of nanoflare amplitudes; 

III. to compare with SUMER transition region and coronal data 
in order to check whether this simple model can reproduce the 
observations and whether the deduced parameters are realistic. 
We begin with a description of the data and their analysis 
(Sect. 2). Next we give an outline of the model used for the 
flare simulation and establish a heuristic reasoning for a log- 
normal distribution of the radiances under the assumption that 
they are entirely due to transient events (Sect. [3]). The results of 
parameter studies with our model are presented in Sect. 4. In 
Section 5 we compare the simulations with SUMER measure- 
ments, and we outline our conclusions in Sect. 6. 



2. The SUMER data: reduction and analysis 

2.1. Data and reduction 

SUMER is a stigmatic normal incidence telescope and spec- 
trometer, operating in the wavelength range from 46.5 to 
161.0 nm, depending on the spectral order and the choice of de- 
tector. For a gener al description of the SUMER instrument and 
its data we refer to lWilhelm etaild 19951) . The SUMER slit with 
angular dimensions of l"x300" used for the data analyzed here 
is imaged by the spectrograph on to the detectors with a reso- 
lution of about 1" per pixel in the spatial direction and 4.4 pm 
per spectral pixel in first order and 2.2 pm per spectral pixel in 
second order. 

The O iv line at 79.0 nm and the Ne vm line at 77.0 nm, which 
are used here, are measured in first order. We also analyzed 



and modelled the Siiv line at 139.3 nm, but the results were 
very similar to those pertaining to the O iv line and are not de- 
scribed in detail in the following. On 8 February 1998, SUMER 
performed a long-duration observation of quiet areas near disk 
centre in the O iv and Ne vm lines. The measurements, which 
have been described by Wilhelm and KalkofenI ( 2003 ). were 
taken over 7 h and 20 min with a cadence of 33.5 s and a par- 
ticularly accurate compensation of the solar rotation. The Si iv 
measurements, recorded 19 July 1998, belong to a quiet Sun 
explosive events study. Recordings were made over 3 h and 
35 min with a cadence of 15 s. (See TableQ]for a summary of 
the measurements.) 

The data were corrected for the Hatfield, the geometric distor- 
tion, and for detector electronics effects such as dead-time and 
local-gain depression. 

After the instrumental corrections and the radiometric calibra- 
tion, the solar radiances were determined by integration over 
the line profiles, which were derived by least-squares fits of sin- 
gle Gaussian functions and a linear background. The spectral 
background (continuum) was subt r acted prior to integration. 
The data used by IPauluhn et al.1 (bOOOh were mainly image 
scans and thus gave the distribution of the radiances in the 
scanned area during the time span of the image recording. They 
may be considered to give a snapshot of the distribution of ra- 
diances. Here, we investigate time series, i.e., we focus on the 
distribution of radiances in a certain area caused by temporal 
variations. 



2.2. Data analysis: lognormal distribution 

Before doing any modelling we investigate whether time se- 
ries of quiet Sun radiances are approximated by a lognormal 
distribution. After finding that this is the case, we determine 
the parameters of the distribution, which is characterised by 
the following probability density function p for an independent 
parameter x (e.g., the radiance) 



p(x) 



No 



;V2lr 



exp(- 



2cr 2 



), 



(1) 



with ju =< log(x) >, cr = -y/Var(log(x)) and No a normalization 
factor. Var means the variance, <x is called the shape parame- 
ter because it determines the shape of the distribution: small <x 
makes the distribution more Gaussian-like, large cr makes the 
distribution more skewed, p. is called the scale parameter and 
stretches the distribution function. 

An example histogram of a SUMER time series of the transi- 
tion region (Siiv line at 139.3 nm) is shown in Fig.Q] along 
with a lognormal fit. The bin size is 5% of the maximum ra- 
diance value. For practically all spatial pixels the distributions 
of the Oiv, the Siiv, and the Nevm radiances are well de- 
scribed by lognormal f unctions. This does n ot automatically 



follow from the result of lPauluhn et alj ([2000), since they con 



sidered the distribution of radiance from different spatial lo- 
cations, which represented a mix of network and cell interior 
locations. Also, the lognormal distributions from network and 
cell deduced from the time series differ in their fi values, but 
not significantly in their <x values (see Sect. [5J. 
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line 


region 


date 


duration/cadence 


remarks 


Oiv 79.0 nm 
Ne viii 77.0 nm 
Si ivl39.3 nm 


TR, T « 10 5 K 
upper TR/lower corona T ~ 7 x 10 5 K 
TR T ~ 10 5 K 


8 February 1998 
8 February 1998 
19 July 1998 


7 h 20 min, 33.5 s 
7 h 20 min, 33.5 s 
3 h 35 min, 15 s 


quiet Sun study, high-telemetry 
quiet Sun study, high-telemetry 
quiet Sun explosive events study 



0.10 



0.08 - 



0.06 - 



0.04 - 



0.02 - 



0.00 




0.0 0.1 0.2 0.3 0.4 0.5 

Fig. 1. Example of a histogram from a SUMER transition re- 
gion time series (quiet area, Siiv line at 139.3 nm). The radi- 
ance data and the corresponding lognormal fit are given by the 
thick-line histogram and diamonds, and a simulation and its fit 
are indicated by the thin-line histogram and a dashed curve. 

At coronal temperatures the shapes of the probability density 
functions are narrower than for the transition region and more 
symmetric (cf., Pauluhn et al., 2000), and the time series shows 
peaks that are less well marked, which will be shown in Sect. [5] 
Thus the difference between the distributions of radiances re- 
flecting different temperatures is equally shown by data sam- 
pling a part of the solar surface or a single point over a time 
interval. The p. and cr values found from fits to the observa- 
tions are therefore useful constraints to any model of quiet-Sun 
radiances. 

3. A simple model 

3.1. Description 

Using a simple model of transient brightenings, we produce 
synthetic time series of EUV radiances. We presume that flar- 
ing is an intrinsically stochastic process, and our radiance vari- 
able is a time-dependent random variable. One simulation thus 
delivers a possible realization of this process. 
Our model consists of a time series of random kicks (acting as 
"flares"), applied to an initial radiance. Each kick is followed 
by the exponential decay of the radiance. The final radiance is 
given by the sum of the radiances of all the overlapping tran- 
sient brightenings. The numerical value of the initial radiance 
is not important, since it is damped just like every brightening 
and the model reaches a statistical "steady state" after a certain 



relaxation time. We have 5 free parameters: the maximum and 
minimum allowed flare amplitude, y max andy mI -„, the power-law 
exponent a, the e-folding or damping time of the flare tj, and 
the frequency of the excitation, i.e., the flaring rate or flaring 
probability p f. This set is complete if we assume that the shape 
of the nanoflare is given by a single kick with a sharp rise and 
a successive exponential decrease. These parameters are suffi- 
cient to quantitatively reproduce the radiance PDFs and power 
spectra. In order to smooth the steep increase a rise time r r 
can be introduced (which, if included, is chosen to be a fixed 
fraction of the damping time and hence does not add an addi- 
tional free parameter). The introduction of a rise time improves 
the qualitative agreement between synthetic and observed time 
series. 

For the energy content of a flaring event of amplitude yo and 
damping time T d the following proportionality holds 

r°° _i^o2 

E ~ I y Q e r " dt = y T d . (2) 
J to 

Thus, the energy content of a flare with amplitude yo, damping 
time Td and rise time r r can be estimated as 

E = qy 7tA(r d + r r ), (3) 

where A is the solar surface area, if we assume equal radiation 
in all directions (Lambertian surface), and q is a fraction of the 
total energy radiated in the observed spectral line. In general 
q is a small number, determining whose exact value is beyond 
the scope of the present paper. 

3.2. The simulation 

The model produces realizations of possible radiance time se- 
ries, i.e., we get the radiance as a time-dependent stochastic 
variable (stochastic process), and it describes an example of a 
process which is close to being Markovian. In a Markov pro- 
cess, the stochastic variable at one timestep t n+1 is only depen- 
dent on the directly preceding timestep t„. It has "nearly no 
memory" of the history, and the probability to find a variable 
at position x n+ i at time t n +\ is calculated by an initial probabil- 
ity and the two-time transition probability (which is the condi- 
tional probability that the variable is in state x„ + \ at f„+i under 
the condition that it has been in state x n at t„). 
The simulation involves the following steps: 

1. Generate a distribution of flare strengths, i.e. (positive) val- 
ues of flare amplitudes f n . 

2. Start from an initial radiance value ro > (the "first kick", 
taken from the flare distribution, ro = /o). 
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3. At random time f, another radiance r, is generated by 
adding a flare kick 



+ fi- 



(4) 



The random amplitude is required to be consistent with 
the chosen distribution function. 
4. At successive time steps tj, j > ;, the radiance values are 



r-j = r r exp( ), 



Td 



(5) 



with Td the damping time, corresponding to a difference 
equation of 



n-n-i = o-i( ex P(- 



Td 



)-l), 



(6) 



multiplicatively generated from the preceding values. Here 
we assume for simplicity that all brightening events have 
the same damping time. 

The probability of a transient brightening occurring (per time 
step), pf with < pf < 1, is simulated by drawing equally dis- 
tributed random numbers between and 1, and a flare event is 
started at f, if the random number falls within a certain fraction 
of the interval (0,1). For example, a flaring probability of 30 % 
or 0.3 is realized by applying the kick if the random number 
is smaller than 0.3. The flare process is thus a Poisson process, 
which is equivalent to the fact that the waiting times, i.e., the 
time intervals between two flares ATf, are exponentially dis- 
tributed with parameter p /, and the mean waiting time < A7" > 
is given by 1/p/- Additionally to the excitation, we include a 
damping process with damping time r^, which is inverse to the 
strength of the damping process. In Fig.|2]an example sequence 
of flares with very large mean waiting times or low frequency 
Pf is shown. For higher frequencies, the radiances due to indi- 
vidual "kicks" or flare-like events overlap much more. 
A heuristic explanation for the distribution of the radiances is 
given in the following. We assume that we can express the kick 
as fi = c - r,_i, with c a positive value, i.e., the new radiance is 



n = r,-_i + c- r,-_i = c r;_i, 



(7) 



a multiple of the predecessor radiance. We stress, however, 
that we prescribe the distribution of the additive components 
f, (given, e.g., by observation as a power law) and not that of 
the multiplicative component c. 
Thus, the radiance value at the nth step is 



ro ■ c\ ■ a 



(8) 



with the c, being random factors if a new kick happens at this 
time. Otherwise c, represents the damping as described by the 
exponential in Eq. ©. 

This assumption of a multiplicative process is plausible in the 
light of the radiance being ultimately generated by emerging 
and decaying magnetic field. This is for example implied by the 
excellent correlation between high-resolution magnetograms 
and images in EUV and X-ray wavelengths, although the un- 
derlying processes are still far from being understood. 




Fig. 2. An example sequence of simulated flares. Here, the flar- 
ing frequency has been chosen to be very small, Pf - 0.08, 
such that the average waiting time between two flares is 12.5 
time steps. 

Taking the logarithm of Eq. ([8]) yields 

log(r„) = log(rn) + log(ci) + log(c 2 ) + • • • + log(c„). (9) 

The Central Limit Theorem then states that as n goes to infinity, 
the distribution of the sum converges to a normal distribution. 
For the Central Limit Theorem to hold, the c, have to be inde- 
pendent random variables, which is not precisely fulfilled here, 
because the flaring process is coupled via the damping time 
scale and thus has a (short) memory. We can, however, assume 
that the damping time is short relative to the other time scales 
of the system, e.g., the duration of the measurements, which 
would link the flare process to a c oloured noise proce ss instead 
of a white noise process (see, e.g jHonerkam p 1990). 
Thus, the distribution of the sum in Eq. (O is approximately 
normal, so that the distribution of the r„ is approximately log- 
normal. 

On the time scales under study (several hours), we can assume 
stationarity of our distribution, which means that (after an ini- 
tial relaxation time) a steady state is reached. This does not 
hold for the solar radiance when the structure of the region 
changes significantly, e.g., due to emergence of new flux and 
thus a change in activity. Our study is thus limited to truly quiet 
regions. 

The two characterizing parameters of the stationary solution 
Eq. {D, p and <x, have been determined empirically by parame- 
ter scans for Gaussian and for power-law stochastic flare input. 

4. Parameter studies 

In this section we describe the effects of varying the free pa- 
rameters of the simulation on the parameters of the radiance 
distribution. The simulated time series consisted of n — 5 x 10 5 
time steps each. 

First we need to demonstrate that the distribution function of 
the radiance has a lognormal shape. We ran the model for flare 
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Fig. 3. Histograms of the simulations with tj =10, 50, 100, 
and 200 time steps for power-law flare distributions with fixed 
flaring probability of 0.2. 



damping is weaker. Then we do the opposite and vary pf 
while leaving unchanged. The results (for the parameters 
a = 2.5,y mi „ = 0.016, y max = 0.8) are shown in Fig. [4] Both 
curves agree with each other within the uncertainties of the fits. 
From Fig. |4]it follows 

V2<x = a + const. (10) 

The proportionality factor K a varies between values of 1 and 2. 
As either t</ or p f increase the radiance enhancements due to 
consecutive kicks increasingly overlap, leading to an increas- 
ingly symmetric distribution. For very low t</ and p f low radi- 
ances dominate the distribution, while the few high peaks make 
it very asymmetric (see Fig. [3]). Also empirically, we find the 
following dependence for p: 

p = \og(T d m f pf) + p f exp(p f ) + const., (11) 
where mj is the mean value of the flare amplitudes. 



amplitudes which are distributed according to a power law, 
with the exponent a — 2.1 and the amplitude ranging between 
0.02 and 3.0 W m _2 sr _1 . Figure|3]shows examples of simulated 
distributions with damping times = 10, 50, 100 and 200 time 
steps, and a flaring frequency of pj = 0.2 (in inverse time 
steps). Clearly, the resulting distributions are highly skewed to- 
wards small radiances for a short damping time. In the limit 
of infinitely short damping time the most common value of the 
radiance is zero. As the damping time increases, so do the av- 
erage value and the width of the radiance distribution, which 
also becomes increasingly symmetric. The same behaviour is 
found as the flaring frequency increases (not plotted). Note that 
for increasing tj, as the overlap between individual brighten- 
ings increases, a background of nearly constant brightness is 
built up. This becomes visible in Fig. [3] through the increasing 
amount of points at small radiances at which the PDF is zero. 
In order to check if a lognormal results only for a power law we 
also considered, as a test case, a Gaussian distribution of flares, 
which results in a far more symmetric distribution of radiances 
than when the distribution of flare amplitudes is a power law. 
Nevertheless, it is also well approximated by a lognormal dis- 
tribution. In the following we consider only the power-law dis- 
tribution of flare amplitudes, since it is supported by observa- 
tion. 

4.1. Variation of the damping time and flaring 
frequency 

The shape of the equilibrium distribution is generally deter- 
mined by the ratio of damping to excitation. Only when the 
two forces balance each other, can a stationary state be reached. 
To confirm this for our simulations, we vary the quotient of 
damping strength and excitation frequency in two dif- 
ferent ways: first, we vary T4 and allow pf to remain con- 
stant. With increasing damping time, the shape parameter cr 
decreases, which means that the skewness of the lognormal de- 
creases. The distribution becomes more symmetric when the 



2.0 



1.5 

V2cr 

1.0 



0.5 



0.0 




0.0 0.2 0.4 0.6 , 0.8 1.0 1.2 

Fig. 4. Variation of the lognormal shape parameter with the ra- 
tio of damping to flaring rate. The stars give the results for con- 
stant flaring rate p / and varying damping time Td, and the dia- 
monds give the results for constant t</ and varying p j. The hor- 
izontal dashed line represents the observed value as deduced 
from SUMER data. 

Where the damping or excitation exceed certain limiting val- 
ues, no stable solution exists that can be described by a lognor- 
mal distribution. If the excitation is too strong (for this set of 
parameters at < 0.1), the realizations grow unboundedly, 

and at the other end of the range, if the product exceeds a 
certain value (here 1 .2), the shape of the radiance distribution 
changes from a two-sided function to a function with a max- 
imum value at or very near zero and a rightwards extending 
tail. 

4.2. Dependence on the minimum and maximum flare 
amplitude 

For constant flaring frequency and damping time the energy in- 
put to the system is determined by the flare distribution, i.e., 
the minimum and maximum flare amplitude and the exponent 
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of the power law. For fixed upper and lower boundaries, the 
energy input decreases with increasing exponent. Thus, in or- 
der to reproduce a certain (measured) radiance time series and 
its mean value, the boundary y m ,„ and/or y max has to be varied 
simultaneously with the exponent such as to keep the energy 
input constant. The mean of the flare radiance input is given by 

m f = E J EE~ a dE = E J E^dE (12) 
with E = E~ a dE^ 

_ (1 GO (ymax ymin ' 

Figure [5] shows the dependence of the lognormal parameters 
on the minimum and the maximum flare amplitude for a 
1 .8. We find that the distribution of the radiances depends more 
strongly on the lower energy limit of the nanoflares than on 
the upper limit also for other a values, including a > 2. The 
mean of the resulting radiance values varies as the mean of a 
lognormal function 

meani ogn = ef + ~. (13) 



computed over the ensemble of all network time series in the 
observed quiet-Sun area (60 pixels). They have been calcu- 
lated with a binsize of 0.01 WrrT 2 sr _1 . The exponent a = 2.5 
yielded the best match to the selected SUMER data. The qual- 
ity of the fit was determined from a comparison of the density 
functions fitted to the histograms and the corresponding cumu- 
lative distribution functions. However, this diagnostic does not 
allow values of a > 2.5 to be distinguished from each other un- 
less the quality of the data is high and the statistics are good. In 
Fig.|7]we show the lognormal parameters p and cr as a function 
of a for constant flare frequency pf = 0.43 (1.3 x 10~ 2 Hz), 
damping time r,/ + r r = 5.72 time steps (191.7 s) and constant 
upper amplitude value of y max = 0.8, the value of a has been 
varied while also varying the y m ;„ value (y,„i„ = 0.002, 0.006, 
0.011, 0.016, 0.019, 0.021, respectively), such as to keep the 
energy input constant (as well as the number of flares). Clearly, 
the distribution of radiances becomes more symmetric when 
a increases (when the amplitude range is selected in order to 
keep the mean value constant). It is also evident from Fig. [6] 
that the sensitivity of the radiance distribution to variations of 
a is highest for small a values. Figure[7]confirms the result ob- 
tained from Fig. [6] that a — 2.5 agrees best with the SUMER 
data. 







a = 1.8 


i 


(a) 





0.000 0.002 0.004 0.006 0.008 0.0 0.2 




(c) 



0.000 0.002 0.004 0.006 




Fig. 5. Variation of the lognormal scale and shape parameters 
with the lower (a,c) and upper (b,d) cutoffs of the flare ampli- 
tude distribution. The values employed for the other parameters 
are a = 1.8, = 8.01, p/ = 0.60, (y„„>, variable, y max = 0.5); 
and (y max variable, y,„j„ = 0.002). The horizontal dashed lines 
represent observed values (averages over Oiv quiet network 
area). 



4.3. Dependence on the power-law exponent 

Figure [6] shows the probability densities and cumulative distri- 
bution functions (CDFs) of the simulated time series for differ- 
ent values of a. Also plotted are the same quantities from the 
time series of the transition region O iv 79.0 nm line recorded 
by SUMER in the quiet network. The histograms have been 



4.4. The power spectrum as an additional diagnostic 

Not surprisingly, just the probability density function alone 
does not lead to a unique solution, since the PDF contains no 
information on time scales. This is illustrated in Figs. [8] [9] and 
[TOl In Fig. [8] we show the observed and simulated PDFs. The 
observed PDF of this network time series is very similar to that 
of the mean which is shown in Fig|6] Obviously, the chosen pa- 
rameters of this particular simulation (a = 2.24, + r r = 19.1 
time steps, pf = 0.16, y m ,„ = 0.014, y max - 0.76) give an excel- 
lent fit to the observed PDF. A comparison with an actual time 
series (Fig. |9]l shows, however, that the individual microflares 
last much longer in the simulations than in the data. Therefore, 
the probability distribution cannot constrain the dynamics of 
the time series. In terms of the free parameters of the simula- 
tion, this implies that multiple combinations of free parameters 
can produce the same probability distribution, i.e., the same p 
and cr values. To demonstrate this we have varied a and — — 

TdPf 

together, such that for smaller a we have chosen a higher — — 
in order to keep the cr near the observed value. The observed 
values of the p and <x parameters are generated using all the 
sets of simulation parameters given in Tab. [2] 
The dynamics of the time series corresponding to the param- 
eters given in Tab. [2] however, are in general not identical. In 
particular, Td and pj influence the duration of the individual 
events and the interval between them. In order to constrain the 
time dependence, a diagnostic sensitive to time scales of ra- 
diance variation is required. The (Fourier or wavelet) power 
spectrum is such a diagnostic. The power spectral density was 
computed here using (Morlet-)wavelets (see, e.g., Torrence & 
Compo, 1998). This choice is driven by the smaller amount 
of fluctuations in this quantity than in Fourier power spectra, 
since wavelets introduce a smoothing in the spectral domain. In 
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Table 2. Values of a, pf, , y mI „ and y max , as well as the p, cr, 
and mean of the distributions produced by them. The values of 
the SUMER time series used for comparison (Oiv network) 
were Hsum = _ 2. 16, cr SUM = 0.98, and mean^y^ = 0.15. 
ACDF gives the summed difference of the cumulative distribu- 
tion functions. 



a 


1.8 


2.1 


2.5 


2.8 


3.0 


Tj+T r 


8.01 


6.87 


5.72 


5.66 


5.55 


Pf 


0.60 


0.52 


0.43 


0.43 


0.42 


1 


0.21 


0.28 


0.41 


0.41 


0.43 


T dl'f 












y miii 


0.002 


0.006 


0.016 


0.019 


0.021 


ymax 


0.5 


0.6 


0.8 


0.8 


0.8 


P 


-2.17 


-2.21 


-2.12 


-2.12 


-2.16 


<r 


0.95 


0.97 


1.0 


0.97 


0.95 


mean 


0.15 


0.14 


0.15 


0.15 


0.14 


ACDF 


0.007 


0.008 


0.004 


0.007 


0.013 




0.5 ; 
0.0 . 



Fig. 7. Variation of the lognormal parameters with the expo- 
nent a of the flaring input distribution. The total energy input 
has been held constant by varying the lower boundary of the 
flare amplitudes accordingly. The dashed lines show the cor- 
responding lognormal-fit values to the averaged SUMER Oiv 
network data. 

Fig. [10] we note that the main peak of the global wavelet power 
spectrum of the observed time series lies at nearly one half of 
the frequency of that of the model. Thus, this diagnostic nicely 
complements the radiance distribution. As we show in Sect. [5] 
the combination of distribution function and power spectrum 
leads to much stronger constraints on the free parameters of 
the model. We cannot rule out, however, that the distinct peaks 
seen in the power spectra are partly due to the restricted length 
of the time series and the thereby limited statistics. Moreover, 
some properties of the flaring process are best seen by directly 
comparing with an observed time series, although this compar- 
ison remains qualitative. 

5. Application of the model to SUMER transition 
region and coronal data 

SUMER data of the EUV lines O iv, Si iv and Ne vm were used 
to constrain the free parameters of the model. These lines repre- 
sent temperatures of the lower (cooler) transition region (O iv, 
Si iv) and the upper (hotter) transition region or lower corona 
(Ne vm). We restricted ourselves to areas where no significant 
change of the activity was discernible, i.e., truly quiet regions 
with no detectable trend in the average radiance during the 
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Fig. 8. The histograms corresponding to the time series of 
Fig-Hi The simulated histogram is represented by the thin solid 
line. The dashed line indicates the lognormal-based fit to this 
distribution. The SUMER histogram is given by the thick line 
and its fit is indicated by the stars. 
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Fig. 9. Example of a simulation featuring a similar lognormal 
distribution as a SUMER (Oiv network) time series but a very 
different temporal variation, so that the time series look rather 
different. Parameters of the simulation are a = 2.24, Td + r r = 
19.1 time steps, pj = 0.16, y m ,„ = 0.014, y max = 0.76. The 
corresponding histograms are shown in Fig. [8] 



measurements. The parameter values obtained from these data 
are listed in Tab. [3] This table summarizes the optimum pa- 
rameters for our model to represent the averages over the cell 
and network time series as suggested by comparison with the 
SUMER data. The summed difference of the cumulative distri- 
bution functions (compare Fig. |6]l is given in the last column. 
The model parameters were varied until a good (statistical) 
match was obtained to the SUMER time series, its radiance 
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0.0000 0.0005 0.0010 0.0015 0.0020 0.0025 0.0030 

freq / s" 

Fig. 10. The global wavelet power spectra of the modelled 
(solid line) and the observed (dashed) time series correspond- 
ing to Figs. [9] and [8] 



distribution and its global wavelet power spectrum. We used a 
statistical power-law distributed flaring process and introduced 
a possible rise time of a flare r r , in order to better match the ob- 
served time series. The flare brightness was chosen to both, rise 
and fall exponentially. In all cases, 77 was chosen to be 75 % of 
the damping time tj, so that the flaring time scale or "effective 
damping time" was r r + tj. Note that the need for a rise time 
T r comes neither from the radiance distribution functions nor 
from the power spectrum, but rather is indicated by the direct 
comparison between the time series (such as those shown in 
Figs. [T21 to [171). Each simulation consisted of n — 24000 time 
steps. For the radiance histograms, a bin size of approximately 
5 % of the respective maximum value has been used. 
We found that adding Poisson photon noise ( y/counts) to the 
model realizations gives slightly better matches when compar- 
ing the single time series (see Figs. Q~2] and Q~3] and[T6land[T7li 
and their histograms. 

5.1. Transition region 

Figure QT| shows two sample histograms from one-pixel time 
series (i.e., from a l"xl" area on the solar surface) of the 
SUMER Oiv data, the first from a darker area of the images, 
corresponding to the interior of a supergranule cell, the sec- 
ond from a brighter pixel (network). Also plotted in Fig. [TT] are 
the histograms resulting from simulations with the parameters 
given in the top two rows of Tab. [3] 

The agreement between measured and modelled histograms is 
reasonable. The majority of the differences can be ascribed to 
the better statistics of the model. 

Two examples of SUMER one-pixel time series of the O iv line 
are shown in Figs.[T2]and[T3] together with the corresponding 
simulations. The photon noise was on average 7 % of the signal 
for the O iv data. 

Figure [14] gives the measured and simulated power spectral 
densities, for both, the darker and the brighter region. The mod- 
elled and the observed power spectra now both display peaks 
at similar frequencies. It is noteworthy that the same model pa- 
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Fig. 11. Histograms and lognormal fits of the simulated (thin- 
line histogram and dashed curves, respectively) and SUMER 
(Oiv line at 79.0 nm) measured (thick-line histogram and dia- 
monds) radiances, (a): from a time series in a very quiet area 
(cell interior), (b): from a network sample of the same size 
(l"xl"). 



rameters reproduce both the network and the cell interior data. 
Only the amplitude of the input process had to be adjusted by 
setting the energy of the smallest flares. The large a value of 
2.5 is obtained because the model attempts to reproduce all the 
diagnostics simultaneously and is not restricted to reproducing 
the brightenings but also the background. 

5.2. Corona 

The coronal time series and the corresponding histograms as 
well as power spectra recorded in cell interiors could be well 
represented by a weakly damped and low-frequent flare excita- 
tion, with the parameters given in the third row of Tab. [3] The 
correspondence between model and data can be judged from 
Figs. [T5R [T6l and IT8h. Roughly 9 % photon noise was super- 
posed on the modelled time series on average. 
Satisfactory results for brighter regions could also be obtained, 
as illustrated in Figs.[T5b,[T7land[T8b. The parameters are listed 
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Fig. 12. SUMER time series of a darker ("cell") region (1 pixel) 
in the transition region O iv 79.0 nm line (upper frame), and a 
corresponding simulation (lower frame). 




Simulated liqhtcurve 




Fig. 13. SUMER time series of a brighter ("network") region 
(1 pixel) in the transition region O iv 79.0 nm line, and a corre- 
sponding simulation. 



in the last line of Tab. [3] + r r needed to describe the Ne vm 
time series and power spectra (Fig.[T8l) in the brighter parts of 
the quiet Sun are roughly 2.5 times larger than in the darker 
parts. The product (tj + T r )/pf, however, differs by a factor of 
only 1.5. The difference in this product is sufficient to explain 




0.0005 0.0010 0.0015 0.0020 0.0025 0.0030 

freq / s" 1 
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0.0010 0.0015 
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Fig. 14. Wavelet global power spectra for the time series shown 
in Fig. [12] (a), and Fig. [13] (b). The solid lines show the simu- 
lated spectrum, and the dashed lines give the spectrum of the 
SUMER time series. 

the relatively small difference in the mean intensity (Tab. 01, so 
that the same y m j„ is used for bright and dark regions. 

5.3. Physical parameters 

Clearly, our model is too simple to give, on its own, an accu- 
rate energy range of the flaring events that reproduce the data. 
However, with some restrictive assumptions a limited estimate 
can be made using Eq. ©. 

The amplitudes listed in Tab.[3]correspond to flare energies be- 
tween 6.4 x 10 16 erg and 1.3 x 10 19 erg for Oiv and between 
4.8 x 10 16 erg and 1.6 x 10 17 erg for Ne vm for events of 1 s 
damping time covering an area of 1 square arcsecond (i.e., en- 
ergy fluxes in the range between 12.6 and 2520 erg s _1 cm -2 for 
O iv and between 9.4 and 31.5 erg r'crn" 2 for Ne vm). 
One grossly simplifying assumption that has already entered 
Eq. © concerns the geometry (radiation from a flat surface 
with no centre-to-limb variation). Another is that all events 
cover the same solar surface area, A. A third is that all events 
have the same temperature (or rather the same temperature- 
density relationship), so that a single q value can be employed 
for all. We do not attempt to compute q, so that we only esti- 
mate the energy radiated in the observed spectral line. For the 
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Table 3. Best-fit parameters for the SUMER measurements of 8 February 1998. 
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T d + T r 
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ymin 
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mean 
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mean s um 


ACDF 


dark px. 79.0 nm 


2.5 


5.721 (191.68 s) 


0.43 


0.004 


0.80 


-3.54 


1.08 


0.042 


-3.58 


1.09 


0.039 


0.0045 


bright px. 79.0 nm 


2.5 


5.721 (191.68 s) 


0.43 


0.016 


0.80 


-2.12 


1.0 


0.15 


-2.16 


0.98 


0.15 


0.0043 


dark px. 77.0 nm 


2.5 


44.40 (1487.5 s) 


0.11 


0.003 


0.01 


-3.67 


0.47 


0.027 


-3.67 


0.45 


0.028 


0.0005 


bright px. 77.0 nm 


2.5 


18.28 (612.5 s) 


0.43 


0.003 


0.01 


-2.86 


0.39 


0.06 


-2.84 


0.39 
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0.0022 
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Fig. 15. Histograms and lognormal fits of the simulated (thin- 
line histogram and dashed curves, respectively) and SUMER 
(Ne vm line at 77.0 nm) measured (thick-line histogram and 
diamonds) radiances, (a): from a time series in a dark pixel 
(l"xl"), (b): from a bright pixel. 



area A we assume a size of 10 13 m 2 , w h ich w ould roughly be 
20 SUMER pixels dAschwanden et alJ d200d) found sizes of 
4 - 200 x 10 12 m 2 using TRACE data), and also adopt durations 
of 200 to 2000 seconds. Under these assumptions we obtain a 
total range of energies between 2.52 x 10 20 erg to 5 x 10 24 erg 
for individual events in the Oiv line and 1.88 x 10 20 erg to 
6.3 x 10 21 erg for events in the Nevm line. Of course, this 
energy radiated away by only single spectral lines is still too 
low by several orders of magnitude to make up for the up- 
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Fig. 16. SUMER time series of a darker region (1 pixel) in the 
upper transition region/coronal line 77.0 nm (Nevm ), and a 
corresponding simulation. 



per chromospheric or co ronal energy losses given, e.g., by 



Withb roe & Noy es ( 1977). 



From the flaring frequencies given in Tab. [3] again assuming 
an event size of 20 SUMER pixels, one can obtain a rough es- 
timate of the number of such events on the total solar surface. 
To reproduce the SUMER time series in the transition region 
line and bright corona including the background, it would need 
roughly 7600 events per second. 2000 events per second can 
generate the radiation in the darker areas of the Nevm line. 
Given all the uncertainties entering the above estimates we pre- 
fer to restrict ourselves to order of magnitude estimates: our 
model suggests a nano-flaring rate of 10 3 to 10 4 events per sec- 
ond. These values can be compared with the number of blinkers 
estimated for the total solar surface (50 events/s, Brkovic et al., 
2000) and that of explosive events (500 events/s, Innes et al., 
1997). 

To estimate the total available flare input for heating one would 
need to know or extrapolate these values to the entire spec- 
trum, assuming as precise areas as possible (geometry, filling 
factor) and prescribe which fraction of the total released energy 
is radiated away (and not, e.g., employed to accelerate the solar 
wind). Here we take a simpler approach. In order to get an idea 
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Fig. 17. SUMER time series of a brighter region (1 pixel) in 
the upper transition region/coronal line 77.0 nm (Ne vm ), and 
a corresponding simulation. 



of the energy released by a typical brightening, we use the to- 
tal coronal radiative energy loss of the quiet Sun F r = 10 5 erg 
1 (corresponding to 6.0786 1 X 10 27 erg s - 1 over the 



cm 

entire solar surface area) given by IWithbroe & Noyesl (119771) . 
Assuming 7600 events per second, we get a mean radiative loss 
of 8 x 10 23 erg per event which is an integral over the full spec- 
trum and an average over the darker and brighter network areas. 
If we further assume the same power-law exponent a = 2.5 
over the entire spectrum, we can estimate the lower energy lim- 
its of the flaring events obtained from our analysis from Eq.[T3l 



m aU 



,0', 



-0.5 

max 



.3,-0.5) 
nun ' 



(y, 



1.5 _ -1.5 

max y m i n 



8 X 10 23 ere 



(14) 



for 7600 events per second. Here we have further assumed that 
ymax » ymin- For 2000 events per second (as was estimated for 
the darker parts of the corona with smaller flaring frequency 
and much larger damping time), the average flare energy would 
approximate 3.04 x 10 24 erg. These values would make the 
brightening events in our analysis nano- and microflares, with 
the smallest ones being rather small nanoflares or large pi- 
coflar es, depending on where the lower energy limit for flares 



Parker dl988l) estimated values of around 1.6 x 10 24 erg 



is set. 



per nanoflare event. 



6. Discussion and Conclusions 

In this paper we consider the diagnostic content of radiance 
time series and simulate them with the help of a simple stochas- 
tic model. The radiances of time series measured by SUMER in 
the quiet Sun are distributed following a lognormal function for 
every spatial point. The properties of these distributions reflect 
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Fig. 18. Wavelet global power spectra for the time series shown 
in Fig. [16] (a), and Fig. [17] (b). The solid lines show the simu- 
lated spectrum, and the dashed lines give the spectrum of the 
SUMER time series. 



those found by IPauluhn et al. ( 2000h in the sense that transi- 
tion region radiation shows broader, more asymmetric distri- 
butions, while the coronal radiance exhibits a narrower, more 
symmetric distribution. These results show that quiet-Sun radi- 
ances exhibit a logno rmal distribution irre spective of whether 
we sample spatially (Paulu hn et al. 2000h or temporally (this 
paper), although in the former case we are mixing radiation 
from network and supergranule cells, while in the latter case 
we are considering these features separately. 
The shapes of the lognormals (Figs.fTT1and[T5]) emphasize the 
lack of points around a radiance of zero, with lines at higher 
temperatures showing this effect more strongly. In the past such 
a behaviour has often been interpreted to indicate the superpo- 
sition of a more or less steady background and overlying peaks 
produced by discrete energy releases, e.g., by microflares. Here 
we assume that all the radiance in the observed lines is caused 
by micro- and nanoflares, so that any background is due to a 
superposition of many small brightenings. 
We have developed a model which attempts to explain the prop- 
erties of the observed time series assuming that all the radiation 
coming from the transition region and the corona is due to the 
cooling of gas that had been heated by nanoflares. From the 
multiplicative version of the central limit theorem it follows 
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that the product of many independent, identically distributed, 
positive random variables has approximately a lognormal dis- 
tribution. (For a discussion of the special features of the lognor- 
mal distribution function see also, e.g., Limpert et al., 2001.) 
For a lognormal distribution to form in a simulated time series, 
the input process has to be positive definite (thus, a Gaussian 
with zero mean as input distribution does not produce a log- 
normal). By the positive input process, the asymmetric shape 
is ensured. 

The shape and scale of the resulting lognormal are strongly 
dependent on the ratio of the damping to the excitation pro- 
cess. The higher the flaring frequency and the longer the damp- 
ing time, the more symmetric the resulting lognormal becomes 
(i.e., the cr values which determine the shape of the lognormal 
become smaller and the distribution becomes more Gaussian- 
like) and vice versa. The scaling of the distribution, which is 
its extension in (peak-) height and width and is given by the 
parameter fi of the lognormal, increases with damping time 
(i.e., with decreasing damping) and with flaring frequency. 
Consequently, if many small brightenings overlap, they pro- 
duce a stronger "background" emission than fewer brighten- 
ings lying further apart. It has to be noted that the SUMER 
time series could be reproduced in a model employing the same 
damping time for weak and strong flares. 
An important part of this work has been to identify represen- 
tations of the data which can serve as diagnostics to constrain 
a microflare model and to determine the free parameters, in 
particular the power-law exponent a of the flare amplitudes or 
energies. The combination of the power spectrum of the ra- 
diance time series and the distribution function, together with 
the time series itself are found to constrain the free parameters 
of the model, including a, relatively well. The most important 
result of this paper is that simulations that best reproduce the 
data had power-law exponents a of around 2.5. We stress that 
we have not carried out a thorough uniqueness analysis here, 
although numerous attempts to reproduce the data with a < 2 
were unsuccessful. Note also that this value is obtained assum- 
ing that all (nano-)flares have the same geometry and cover 
the same area. The large a is determined by the need to re- 
produce the background emission seen in the transition region 
and the corona. Smaller a values give too much weight to single 
larger flares which do not reproduce the purely quiet time series 
(cf., Fig. |9). By assuming that all micro- a nd nanoflares cover 



the sam e surface area of 10 13 m 2 following lAschwanden et al 
feOOOh . we find a nanoflaring rate of 10 3 - lOV 1 . 
An interesting result of the current investigations is that the 
damping time scales of the simulations reproducing the coro- 
nal time series are found to be significantly larger than in the 
transition region (by a factor of 4 to 8). A simple estimate of 
the radiative cooling times based on the Rosner formula gives 
(cf., Rosner et al. 1978): 



Trad = 3 



1.38 x 10" 



(15) 



For the temperature of formation of the O iv line, T = 1 .5 x 



10 5 K, £ = 6.3 x 10 



-22 



= 0, while for T = 6 x 10 5 K, 



appropriate for Ne vm, £ = 3.98 x 10 n , a — —2. If we assume 
a constant pressure of 0.3 dyn cirT 2 (= 0.03 N irT 2 ), we obtain 



a ratio of r ra j(Ne vm)/T ra j(0 iv) « 100, which is in the same 
direction, as that obtained from the analysis, although an order 
of magnitude larger. 

The similar flaring probability, pj, obtained for (the brighter 
part of the quiet) corona and transition region does not allow us 
to determine whether the brightnings in the corona and the tran- 
sition region are due to the same or to different flaring events 
based on this analysis. 

Further work is needed, both to demonstrate the uniqueness of 
the deduced parameters, in particular a, to increase the amount 
of realism in the model, and to consider more data. The first 
may be achieved by automating the process of finding a best 
statistical fit to the data (i.e., best simultaneous fit to the dis- 
tribution function and the power spectrum), possibly coupled 
to a Monte Carlo simulation or involving genetic algorithms. 
An analysis of lines sampling higher temperatures is required 
in order to reveal the behaviour of hot coronal gas. There is 
also a need to statistically simulate time series of images. Such 
data may also provide some information on the size distribu- 
tion of flares. Since a single time series belonging to a par- 
ticular pixel may sample only some spatially limited parts of 
particular brightenings, the parameters obtained by statistically 
reproducing it with a 1-d model need not be entirely represen- 
tative of the real Sun. Thus, higher-dimensional models could 
be of great advantage. 

Finally, an application of the analysis developed here to active 
regions could reveal similarities and differences between the 
flaring processes in active regions and the quiet Sun. 
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Pauluhn & Solanki: A nanoflare model of quiet Sun EUV emission 




Fig. 6. Variation of the radiance distribution (the radiance probability density function, shown on the left, and its cumulative 
distribution function, shown on the right) with the exponent a of the flare distribution. The following parameters have been kept 
constant: pf = 0.43, Td + r r - 5.72 time steps, and y max = 0.8. The simulated histogram is represented by the thin solid line. The 
dashed line indicates the lognormal-based fit to this distribution. Also shown are a SUMER O iv histogram representative of the 
network and the corresponding lognormal fit (stars). 



